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SUMMARY 


This paper develops a near-optimal guidance law for generating minimum fuel, time, or cost fixed- 
range trajectories for supersonic transport aircraft. The approach uses a choice of new state variables 
along with singular perturbation techniques to time-scale decouple the dynamic equations into multiple 
equations of single order (second order for the fast dynamics). Application of the maximum principle to 
each of the decoupled equations, as opposed to application to the original coupled equations, avoids the 
two point boundary value problem and transforms the problem from one of a functional optimization to 
one of multiple function optimizations. It is shown that such an approach produces well known aircraft 
performance results such as minimizing the Brequet factor for minimum fuel consumption and the energy 
climb path. Furthermore, the new state variables produce a consistent calculation of flight path angle 
along the trajectory, eliminating one of the deficiencies in the traditional energy state approximation. 
In addition, jumps in the energy climb path are smoothed out by integration of the original dynamic 
equations at constant load factor. Numerical results performed for a supersonic transport design show 
that a pushover dive followed by a pullout at nominal load factors are sufficient maneuvers to smooth 
the jump. 


INTRODUCTION 

The purpose of this work is to develop and implement a near-optimal guidance law for use in 
an aircraft synthesis computer code, such as the ACSYNT code 1 developed at NASA Ames Research 
Center. Of primary interest is the optimization of supersonic transport trajectories. ACSYNT, like 
other such codes, models all aspects (aerodynamics, propulsion, structures, weights, etc.) of an aircraft 
design to produce consistent performance estimates. It is capable of computing “closed” vehicles, that 
is, designs that meet mission requirements, by iteratively adjusting vehicle parameters. It is also capable 
of optimizing design parameters, again by iteratively cycling through the code. 

A key element of any vehicle synthesis code is the trajectory calculation. Because the trajectory 
routine is exercised repeatedly in the course of a design study, it must be efficient, robust, and user- 
friendly. “Exact” trajectory optimization, relying on optimal control theory, requires iterative solution of 
an unstable two-point boundary-value problem (2PBVP), and therefore is not suitable for this application. 
Thus simplifying approximations are required. 

It has long been known that if there is but one state equation, then the functional optimization 
problem (2PBVP) reduces to a function one^’*^. A natural and well-established way to effect the 
required order reduction is to time-scale the system state equations and then apply singular perturbation 
techniques (see for example refs. 2-9). If each state variable is put on its own time-scale then the 
problem is thereby reduced to a sequence of function optimizations. 

The main problem with completely time-scaling the aircraft dynamics is that speed and altitude 
are not time-scale separable. This is usually resolved by replacing the speed by the total mechanical 
energy as a state variable (see for example refs. 2-12), and we adopt this approach here. In addition, 
another new state variable is introduced to replace the altitude, one which removes the inconsistency in 
flight path angle 11,1 ' 1,14 that occurs in the energy dynamics with the usual formulation. This does not 



directly impact the energy dynamics solution but increases the accuracy of the altitude/fiight path angle 
dynamics solution. 

The energy-state approximation (neglecting all dynamics except the energy dynamics) has been 
applied with success to a wide variety of aircraft, including high performance supersonic aircraft and 
launch vehicles. It is perhaps best suited, however, to transport aircraft because the benign maneuvers of 
these vehicles make the assumptions involved in the energy-state approximation (ESA) less questionable. 
The ESA has been applied most thoroughly to subsonic transport aircraft by Erzberger 15-1 . The results 
were so satisfactory that the resulting algorithms currently are being used for on-board guidance in 
commercial transports. 

Applying the ESA to supersonic transports introduces some new features. First, these aircraft have 
higher speeds and usually longer ranges than do subsonic aircraft. More importantly, due to the rise in 
drag near transonic speeds, they typically have an instantaneous altitude change in their energy-climb 
paths. These altitude jumps have been investigated by various means in references 12, 18-22. In this 
paper, we use the approach of references 12 and 22 to address this problem. 

Finally, some numerical results are presented to demonstrate the utility of the method. 


DYNAMIC MODELING 


Equations of Motion 

The equations of vehicle motion in ACSYNT are: 

m = —(5 — —nT m C 

x = v cos 7 

T cos a — D — mg sin 7 
v 

m 

h = v sin 7 

T sin a + L — mg cos 7 

7 ~ 

mv 


( 1 ) 


These equations assume no winds, thrust direction fixed with respect to the aircraft body, and a 
non-rotating flat earth. A linear throttle is not assumed; that is, specific fuel consumption, C, varies 
with thrust. The symbols used here and throughout this report are defined in Appendix A. 


To simplify the terms, define the tangential and normal load factors as 

(Tcosa — D) 

~ rng 


N = 


(T sin a - 1 - L) 

mg 


( 2 ) 
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Then equations ( 1 ) become 


m = —j 3 
x = v cos 7 

v = g(F- sin 7) (3) 


h = v sin 7 

7 = -(AT — cos 7) 
v 


In ACSYNT, as in many other vehicle synthesis codes, equations ( 3 ) are numerically integrated, 
with the 7 term set to zero, for a specified set of ordered pairs of altitude and Mach number (or speed). 
The methods used for this integration are given in Appendix B. The ( h,M ) points needed for the 
integration may come from any number of sources, for example a constant dynamic pressure (constant 
equivalent airspeed) path or an external trajectory optimization. It is our purpose to develop an algorithm 
that generates these points near-optimally for some prescribed cost functional “on the fly”, that is, as 
the trajectory integration proceeds. 


Transformation to New State Variables 

Experience has shown that the state variables in equations ( 3 ) have a natural time-scale separation 
for most vehicles and most missions, except that h and v are on almost the same time scale. To time 
scale separate these variables, we seek a new variable, E(h , v), to replace v, such that the state equation 
for E is independent of 7 11 - 1314 . Taking the time derivative of E and using equations ( 3 ): 

E = E it h + E v v = Eh v sin 7 + E v g(F - sin 7) ( 4 ) 


Throughout, the following notation will be used: If Q is any function of h and v, then 


Qh 



Q = 9Q 

Qv dv 


h 


( 5 ) 


If E is to be independent of 7, from equation ( 4 ): 

E h V - Ey g = 0 


The solution of this equation is 


E = h + 



(6) 
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or any once-differentiable function of this. From equation ( 6 ) we see that E is just the total mechanical 
energy of the aircraft per unit weight. Substituting E for v as a state variable gives 


771 = —(3 
x = v cos 7 

E = vF ( 7 ) 

h = vsin7 

7 = -(N - COS 7 ) 

v 

Numerous analyses have shown that there is a strong time-scale separation between E and h (see for 
example refs. 23, 24). In equations ( 8 ), v is to be regarded as a function of E and h, as given by 
equations (7): 

v = j2g(E - h) ( 8 ) 

The product vF in the third of equations (7) is usually called the specific excess power. 

Equations (7), along with a suitably defined cost functional, define an optimal control problem in 
the five states m, x, E, h, and 7 , with control a (and possibly throttle if it is allowed to vary). The 
boundary conditions on these states are 

m(0) = 7770 

x(0) = 0 
E{ 0 ) = E 0 
MO) = h 0 

7(0) = 70 

where tf is free. 

The following constraints are placed on the trajectory: 

1. Maximum dynamic pressure, q(h,v ) < q m 

2. Maximum Mach number, M(h,v ) < M m 

3. Maximum lift coefficient, ci(h,v) < c^ m 

4. Minimum terrain limit, h> h m 

5. Maximum loft ceiling (locus of flight conditions for which F = 0 for maximum throttle and 
7 = 0) 


m(tf) free 
x(tf) = R 

E(t f ) = E f (9) 

h(tf) = hf 
7 {tf) = If 
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All of these constraints may be written as functions of h and M or of h and E\ when drawn in 
the ( h,M ) plane (fig. 1) they define the flight envelope. In the context of equations (7) they are state 
inequality constraints of the form: 


«t(M)<0; i = 1, 5 (10) 

Optimal control problems with state inequality constraints are a difficult class of problems for several 
reasons 25,26 . 

The complete time-scale decoupling of equations (7) will be formulated later. At present, for the 
sake of dynamic modeling, it is instructive to consider the energy-state approximation (ESA) associated 
with equations (7); it is: 


m = const 
x = const 


E — vF 


( 11 ) 


0 = i;sin7 

0 = -(N - cos 7 ) 
v 

The fourth of these implies that 7 = 0 and the fifth then gives a as a function of h and E. The problem 
thus reduces to a single state equation with h (and possibly throttle) as control and E as state. The 
solution, for a suitable cost functional, may be put into the form (see later) 

f(h,v) = 0 ( 12 ) 

This will be called the energy-climb path, or ECP. This may be either one of the constraints equations (10) 
or an interior extremal. One of the main advantages of the ESA is that it converts the state variable 
inequality constraints, equations ( 10 ), into state-dependent control inequality constraints, a much simpler 
situation from an optimal control point of view. 

Since equation (12) generally gives h ^ 0, 7 will not be zero on the ECP, giving a contradiction. 
What is needed is a new variable that is constant along the ECP. An obvious choice is / itself since by 
equation (12) / = 0 along the ECP 11,13,14 . Since df = dh + f v dv = 0 we have 

— = -llL (13) 

dh f v 


But from equations (3) 


dv _ g(F — sin 7) 
dh vsin 7 


(14) 
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so that 


7 = sin 


l 



This is the consistent value of 7 along the ECP. Also, from equations (3) 

f = fhh + fvV 


(15) 


/ = //i v sin 7 + /u g{F — sin 7 ) (16) 

Note that the choice of variable / actually depends on the nature of the ECP and may vary along the 
trajectory. The equations of motion in the new variables are now: 

m = —(3 
x = v cos 7 

E = vF (17) 

/ = f h vsin-f + f v g(F — sin 7 ) 

^ = ^(jV - COS7) 

V 

These equations are entirely equivalent to equations (1). 

Some examples of the function / will now be given. 

1 . ECP on a terrain limit: 

f(h,v) = h - h m = 0 


From equations (15) and (16): 


fh = 1 , fv = 0 


7 = 0 
f = v sin 7 


2, ECP on a dynamic pressure limit: 

f(h,v) = ^p{h)v 2 - q m = 0 
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fh=2 PhV ’ fv = P v 


7 = sin 


/ = \ph vi sin 7 + pvg(F — sin 7 ) 

3. ECP an interior unbounded extremal. In this case ( d(vF)/dh)£ = 0 so that: 

v 2 

f(h , v) = F + vF v F h = 0 

9 

v 2 

fh — F h "i" vF vh g F hh 
2 u v 2 

fv = 2 Fy + vFyy F fr F fry 

9 9 


-1 




t 


7 = sin 1 


V 1 g(2gF v +gvF vv -2vFf l -v' z F vh J 

( v 2 \ ( 2v v 2 \ 

f = v f F h + vF vh - ~ F hh) sin 7 + g I 2F V + vF vv - —F h - ~ F hvj ( F ~ sin l) 




Since this latter case involves second derivatives of F, usually a severe problem when dealing with 
numerically defined functions as in the case here, for this case it is probably preferable to compute 7 
along the ECP directly from equation (14) 


7 = sin 


1 


F 


1 + 


dv 

Wi 


(18) 


where dv/dh is evaluated numerically along the ESA solution. These examples show that the usual 
choice of variable in the ESA, h, is only valid when the ECP is on a terrain limit. 


Now consider the ESA associated with equations (17): 


m — const. 


x = const. 


E = vF 


0 = fhv sin 7 + f v g(F - sin 7 ) 

0 = -(N - cos 7) 

v 


7 



The fourth and fifth of these are to be solved for ct and 7 as functions of E and /. Direct elimination 
of 7 gives 

0 = fh v (a / 1 ~ N 2 ) + fv9 (-F “ \/l - N 2 j 

and thus the restriction -1 < iV < 1 must be imposed. Since we will need to consider cases N > 1 
later on, this restriction is unacceptable. The problem is resolved by making the small 7 assumption 
(sin 7 = 7 , cos 7 = 1 ), a very good approximation for transport aircraft whose flight path angles are at 
most a few degrees. Finally then, the equations of motion we shall be dealing with are 

m = —(5 


with boundary conditions 


x = v 
E = vF 

f = fhVl + fv9(F ~ 7) 


7 = ^- 1 ) 


(19) 


m(0) = mo 

m(tf) free 

x(0) = 0 

CS- 

II 

& 

II 

o' 

C“+. 

II 

s 

3 

II 

c-+- 

II 

II 

O 

l(tf) = 7/ 


( 20 ) 


where the boundary conditions on E and f are determined by the boundary conditions on h and v, and 
tj is free. 
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OPTIMAL CONTROL AND SINGULAR PERTURBATIONS 


The Maximum Principle 

All of the equations of motion of the previous section (see for example equations (19)) are of state 
variable form: 


x = f[x, u) 


( 21 ) 


where x G JR n is the state and u G U C lR m is the control. Suitable boundary conditions on the state 
vector components are prescribed (see for example eqs. (20)). It is desired to find the components of y 
along the trajectory such that a cost functional 


J = 



( 22 ) 


is minimized. It is assumed that the final time, tj, is free. Extensions of this basic problem such as for 
terminal cost or fixed final time are easily made, but are not of interest here. 


Theorem (the maximum principle) 25,27,28 : Introduce the variational Hamiltonian function 


# = A o 0+f> ifi 


2—1 


(23) 


where the components of the adjoint vector. A, satisfy the differential equations 

; _ dH 

A { — j i> 3. , ■ • , 

Then, if u is an optimal control, there exists a nontrivial solution of equations (24) such that 


(24) 


(a) u — arg max H 

uel 7 

(b) H = 0 

(c) Transversality conditions (“natural” boundary conditions on the A^) hold 

(d) Ao = const. < 0 

In the sequel it is assumed that Ao = 0 does not lead to a solution and therefore we may take 
Aq = — 1 (this scales the adjoint variables A j). 

The maximum principle gives the control as a function of time or of the state variables. When this 
function is substituted into equations (21) and (24), the result is a 2 n dimension 2PBVP in the states 
and adjoints. Exactly n boundary conditions are provided at t — 0 and the other n at t = tf (due to 
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the transversality conditions). Further, the equations are unstable in the sense that if they are linearized 
about a nominal trajectory, one-half of the system matrix eigenvalues will have positive real parts and 
the other negative (unless some are zero). Although many approaches have been developed to solve this 
class of problem, they are all computationally expensive (requiring repetitive solution of the equations), 
non-robust (due to the instability), and not user-friendly (requiring extensive input by experts). Thus 
they are unsuitable for use in a vehicle synthesis code and approximations must be developed for this 
purpose. 


Approximation Techniques 

Our basic approach is to reduce the complexity of the trajectory optimization problem by seeking 
means of reducing the problem to sub-problems of lower order. There are two keys observations in this 
regard. 

First, suppose there is a state variable, say Xj, such that Xj does not appear in the system functions 
/ nor the cost function /q, except for possibly fj, and the final value of xj is unspecified. Then from 
equation (24) and the transversality conditions, the differential equation for the corresponding Xj and 
its boundary condition are 

~ ~~dxj ^ = 0 

The only solution to this linear differential equation for a finite value of dfj/dxj is Xj = 0. Thus, from 
equation (23), we see that the ; th state equation does not influence the optimal control; this equation 
has uncoupled from the problem and may be integrated after the optimal control problem has been 
solved. This is the reason, for example, that the range equation uncouples from the other equations in 
the minimum time-to-climb problem. 

Second, suppose that there is only one state equation ( x is a scalar) and one control variable: 

x = f(x, it) (25) 

with cost functional 

J = J * (p(x, u)dt (26) 

We have then, from equations (23) and (24) 

H = -<f> + Xf 

A dx A dx 

The maximum principle gives, assuming that unbounded optimal control exists, 

H = -4> + Xf = 0 


™ = Jjt + x a -L = o 

du du du 
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Eliminating A from these two equations gives 




( 27 ) 


This may be thought of as an equation for u as a function of x, i.e., a feedback control law. 


Alternatively, a direct approach may be used. Combining equations (25) and (26) gives 

J = [ f ~ f dx 
JO f 

Thus (4>/f) is to be minimized with respect to u holding x fixed. Carrying out this minimization for 
unbounded control results in exactly equation (27). Actually, a stronger result holds for the single state 
case; if u is a bounded control of several components, then the optimal control is given by 2,3 


u = arg min 

ueu 



x=consl 


(28) 


Singular Perturbations and Time Scaling 

We have just seen that if the dynamic system can be approximated by a single state equation, or by 
a series of such equations, then the solution may be obtained by elementary means, without solving the 
2PBVP. Singular perturbation theory provides a framework for accomplishing this, and indeed many of 
the references cited in the Introduction use this approach. 

The extensive literature on the application of singular perturbation theory to optimal control prob- 
lems in general and flight path optimization in particular will only be reviewed briefly here. 

Perturbation methods have a long history of application in applied mathematics. Noteworthy 
examples are viscous fluid flow, nonlinear oscillations, and orbital dynamics. Singular perturbation 
methods were put on a solid mathematical foundation for ordinary differential equations by Tikonov 29 
and Vasileva 30 . Initial applications to control were by O’Malley 31 and Kokotovic 32 . The theory 
concerns differential equations which depend on a parameter in such a way that the solutions as the 
parameter tends to zero do not approach uniformly the solution with the parameter set to zero. 

The regions of nonuniform convergence are modeled by “boundary-layer” equations, a term arising 
in fluid dynamics. Solutions in the outer regions (away from the boundary layers) and the inner regions 
(the boundary layers) are independently determined by expanding all system variables in asymptotic 
power series. These solutions are then “matched” to determine their constants of integration. The final 
step is to combine the solutions to give uniformly valid approximations to the solution of the original 
problem. Thus the procedure is termed the method of matched asymptotic expansion (MAE). 

Experience has shown that for the highly dynamic maneuvers of high performance figh ter/ attack 
type aircraft, carrying out the expansions to first order is required for high accuracy (see refs. 7 and 8 for 
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example). For low performance aircraft, such as commercial transports, however, zero order analysis 
has been found to suffice (refs. 15-17 for example). The exception, for supersonic aircraft, is the rapid 
altitude transition typically occurring at transonic speeds; study of this transition is one of the main 
objectives of this report and will be taken up in detail later. 

In this report, for the most part, we will consider only zero-order approximations and complete 
time-scale decoupling. For this simple case the elaborate procedures of the MAE method are trivial® 
and do not need to be further explained. 

Reference 33 was the first to suggest complete time-scale decoupling and to recognize its advan- 
tages. In this approach, a “small” parameter e is inserted into the equations of motion as follows: 

*0 = fo(s.,u) 


e n x n = fn(x,u) 

or 

t l Xi = fi(x,u) ; i = 0, --,71 (29) 

where now x = (x 0 ,zi, ••,:£„). The maximum principle for the system (29) is the same as before, but 
with (see Theorem 5.1 of ref. 8) 

= ( 3 °) 
i = 0 



i = 0,--,n 


(31) 


The i 1 * 1 dynamics are obtained by the stretching transformation t t = t/eK Substituting and then 
setting e = 0 gives (where now the dot denotes differentiation with respect to t{) 

x 0 = 0 => xq = const. 


Xj_i = 0 i = const. 

±i = fi W 

0 = fi + 1 

0 — fn 
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Thus the variables on a slower time-scale than x* are held constant arid the variables on a faster time- 
scale than Xj have their system functions set to zero. In order to be able to apply the maximum 
principle to this single-state problem, the conditions of Theorem 5.3 of reference 8 must hold. Let 
/ = (/ i+1 , ■■■ ,f n ) and xj = (xj+i, x n ). Then the key condition is that the matrix 


' d _k d J=£ 

dx_f ’ du _ 

have maximum rank evaluated along the solution. 


(33) 


If condition (33) is satisfied, then by the implicit function theorem the equations 0 = /y can be 
solved for n — i of the components of xj and u in terms of the remaining m. After substituting these 
solutions into Xj = /j, the optimal control may be determined directly from equation (28) with /j 
replacing /. Alternatively, the equations fl = /y may be adjoined with ordinary Lagrange multipliers 
to the Hamiltonian function and the maximum principle applied. This latter method has the advantage 
that it provides the values of these multipliers. This is of interest because these multipliers are the slow 
estimates of the adjoint variables associated with the fast states 8 . 


In the following section, transport aircraft guidance laws will be developed using the following 
time-scale dynamic model associated with equations (19): 


m = —e/3 


x = v 


eE = vF (34) 

e 2 / = fh v 7 + fvg{F - 7 ) 
e 2 i = -( N-l) 

V 

Note that with this formulation the mass is constant on all time-scales to zero order. The implications 
of this will be discussed later. 

Note also that the system is not completely time-scale decoupled because / and 7 are on the same 
time-scale. This was the approach adopted by Ardema (with h replacing /) 7,8 . Calise, on the other 
hand, time-scale decoupled h and -y2,3,4,33 w jjj 5 e discussed in more detail later. 

As a cost functional, following Erzberger a weighted sum of flight time and fuel consumption is 
adopted 15-17 . 


j= [ tf {Ki + K 2 (3)dt 

J 0 


(35) 
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Since some elements of transport airplane direct operating cost are time dependent and some are 
fuel consumption dependent, a proper weighting of these two effects by appropriate selection of the 
parameters K\ and K2 will give a close approximation of direct operating cost. 

Finally, note that the system dynamics do not depend on state variable x and that therefore the 
state equation x = v would uncouple from the problem if its terminal condition were not specified. 
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GUIDANCE LAW DEVELOPMENT 


Range Dynamics 

Setting e = 0 in equations (34) gives the range dynamics: 

771 = 0 
X = V 

0 = vF 

0 = fh,vi + fvg{F ~ 7) 

0 = -(N - 1) 

v 

Thus the single state equation with its boundary conditions is 

x = v , x(0) = 0 , x(tf) = R 


subject to 


771 = const 

F = 0 


7 = 0 


N = 1 


(36) 


(37) 


(38) 


The matrix (33) evaluated for conditions (38) is 


vFe 

vFj 

0 

vF a 


IlfFe 

UgFf 

fhV - fv9 

fvgFa 

(39) 

9 V 

—Nf 

V J 

0 

-N a 

V J 



where, if Q is any function of E, f, and a. 


n dQ 

Qe =8E 




Qf — 


dQ 

df 


E,a 


Q -*9. 
Qa ~ da 


(40) 


f-.E 


The rank of matrix (39) depends on the energy dynamics solution, which determines / . For example, 
if the energy dynamics solution is on a terrain limit, then / = h — h j* so that fh = 1 and f v — 0. Thus 
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the matrix (39) becomes 


vFg vFj 0 vF a 

0 0 v 0 

-Ne ~N f 0 9-N q 

. v v J v J 

For the special case of thrust-aligned-with velocity, N = LjW and N may be taken as the control; the 
matrix now becomes, with h replacing /, 

vFe vF h 0 vFpj 

0 0 v 0 

0 0 0 - 

v J 

Clearly this will have maximum rank if either Fe ^ 0 or F^ ^ 0. 

Assuming that matrix (39) has maximum rank, we may apply the maximum principle to the single 
state problem defined in equation (37). Although equation (28) could be used to directly determine the 
optimal control, because the adjoint X x will be needed we proceed by forming the Hamiltonian. Note 
that the constraints (10) are now control constraints and do not need to be adjoined to the Hamiltonian. 

Forming the Hamiltonian (see eqs. (23), (35), and (37)): 

H = -Ki-K 2 0 + X x v (41) 

subject to F = 0, N = 1 and equations (10). Applying the maximum principle gives the optimal control 
as 


h c , E c = arg min 
h,E 


^ Ki+K20}j 


F = 0 
N = 1 
eqs. (11) 


and the value of A x as 


A 


X 


Ki + Krfc 
Vc 


Equation (42) defines the optimal cruise conditions. 


(42) 


(43) 
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There are two interesting special cases. First, if K\ = 1 and K 2 — 0, the problem reduces to 

v c = max(u) (44) 



as expected for minimum time. Second, if JiTi = 0 and K 2 = 1 and the fuel-flow varies linearly with 
throttle near the cruise point, equation (42) is equivalent to 


he, Ec — 


arg max 
h,E 


ML/P) 

C 


(45) 


where C is the thrust specific fuel consumption. That is, the Brequet factor is to be maximized. 

The total range of a transport aircraft is the sum of the ranges covered during the ascent, cruise, and 
descent portions of the flight. In our analysis of the range dynamics, the ascent and descent portions of 
the flight occur on a faster time scale and thus do not appear in the determination of the cruise condition. 

In Erzberger’s analysis of this problem 15-17 he subtracts out the range covered in climb and descent 
in determining the cruise conditions. This is important in short range flight and in fact Erzberger was 
able to get good results for flight ranges short enough to be composed entirely of climb and descent. For 
the long range flights of supersonic transports, of primary interest here, this factor is of less importance. 
In the context of singular perturbation theory, climb and descent range may be expected to appear as 
first order corrections. 


The range dynamics solution assumes constant mass. Variations in mass between take-off and 
cruise when determining the cruise point may be expected to be accounted for by first order corrections, 
not pursued here. 


Energy Dynamics 

Changing the independent variable to t\ — t/e in equations (34) and then setting e = 0 gives (the 
dot will denote differentiation with respect to t\ in this section) 


m = const 
x = const 
E = vF 


F 



N = 1 


The matrix (33) for this case is 

A fh v ~ fv9 fv9F a 

?-N f 0 -N a 

. v J v 


(46) 


(47) 
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where 


F r 

A = j _ j \vf hh hf f v g + vg f hv v f f v + f h v f f v g 

, , fv 9 (F h h f + F v v f ) 

~fvh hf g fhV - f vv Vf g f h v + p ' 


For the case of solution on a terrain limit, / = h — hj- , fh = 1, f v = 0 , f^h = fvv = fhv — 0 so that 
(47) becomes 

0 v 0 

-Nf 0 -N a 

v 1 v 

For the special case of thrust-aligned with velocity vector and N replacing a as control, this reduces to 

0 v 0 

0 0 - 

v 

which is in agreement with Section 6.2 of reference 8, and clearly has maximum rank if v ^ 0. 
Forming the Hamiltonian associated with equations (46): 

H = —K\ - K 2 0 + X x v + X e vF (48) 

The constraints (10) are state-dependent control constraints for this problem. Maximizing H gives 


h — arg max 

h 



E = const 


(49) 


where 


P = 


vF 

K\ + K 2 /3 - X x v 


(50) 


with the value of X E as 

A £ = i (51) 

Note that as h and v approach h c and v c , P becomes infinitely large. The three terms in the denominator 
of P have the following obvious interpretation. In climb, three factors are important: minimizing time 
(K\ term), minimizing fuel consumption (K 2 ), and covering range (— A x t>). 
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For the case of an unbounded local maximum, equation (49) implies 



E 


= 0 


or, in terms of v and h (see eqs. (5)), 


Pv- 



( 52 ) 


For this case, 

t(v,h) = P v --P h 
fh = Pvh — 

fv = Pvv Ph Phv 
9 9 


Substituting into the third of equations (46), 


7 = 


F 


1 + 


9i v Pvh+Ph—9P uv) 


(53) 


As mentioned earlier, it is probably best to avoid computing numerical second derivatives and use 
equation (18): 


7 


i , v dv 

1 + 5 m 


instead for the value of 7 along the energy dynamics solution. 


Fast States Dynamics 

Changing the independent variable to <2 = f/e 2 in equations (34) and then setting e = 0 gives (dot 
denotes differentiation with respect to t<i): 

m — const. 
x — const. 

E = const. (54) 

/ = fhV 7 + fv9(P ~ 7 ) 

7 = ~(N — 1 ) 

V 
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Our main interest in this paper is to use these equations to model the altitude transition that typically 
occurs transonically in the energy dynamics solution for supersonic aircraft. There have been three 
approaches to the solution of equations (54). 

Ardema 7,8 for the case of / = h, left h and 7 on the same time scale and iteratively solved the 
associated 2PBVP. Although this is not the approach that will be used here, the problem is formulated 
in general in Appendix C as a starting point for future investigation. Calise 33,34 time-scaled decoupled 
h and 7 and obtained non-iterative solutions for each. This required adding a penalty term on 7 to the 
cost function and a “constrained matching” technique. 


The approach used in this paper is a non-optimal one that assumes the fast state dynamics occur 
at constant load factor, N. This is motivated by reference 22 which showed that the transonic altitude 
transitions occurring in discontinuous energy dynamics solutions consisted of a push-over followed by 
a push-up (see fig. 2 which is reproduced from ref. 22). Reference 12 modeled this load factor history 
by two constant load factor segments and obtained good results. Using a non-optimal approach to the 
fast dynamics is partially justified by the fact that these altitude transitions take relatively little time and 
consume relatively little fuel. 


One way to approximate the altitude transition is to begin flying a constant minimum load factor 
flight path when a jump is detected and then switch to a constant maximum load factor when the new 
branch is crossed; this is the dotted path in figure 3, from reference 12. This is undesirable for two 
reasons. First, the transition is initiated too late, and second, the transition path overshoots the new 
branch of the energy dynamics solution. In our approach, we use the fast state dynamics to determine 
v, the optimum point for transition through E (see fig. 3). 


Noting that F = 0 because E = 0, consider the last two of equations (54) 

/ = fh v l ~ fv91 

. gK 


(55) 


where K = N - 1 is a known constant. 

Following reference 12, the first of these equations is divided by the second to give an equation in 
/ and 7: 


d'y 


~ fvd)l 
9 k 


9K ! ^m = I ldl 


From equation (9) for E = const., dh = — — dv so that 


df = fhdh + f v dv — if v — fh~J dv 


( 56 ) 
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Substituting into equation (56) and carrying out the integration gives 


—K 


/v = / 7<i7 


+ const. 


1 9 

—K In v — -7 + const. 
£ 


( 57 ) 


Now label the last point on the subsonic climb path as point 1, the first point on the supersonic 
climb path as point 2, and the load factor transition point by an overbar (see fig. 3). Then equation (57) 
must hold from point 1 to the transition point with K\ — N\ — 1 and from the transition to point 2 with 
K 2 = N 2 -1: 




Solving for v and 7: 



(58) 


7 A 


\ 


2 (K^l - K 2 ~l) + 4 A'jA'2 In g 

2(/fi - K 2 ) 


This is the same solution as obtained in reference 12 except that now the values of 71 and 72 are to be 
determined according to equation (15). 

The transition path is then determined as follows. Constant load factor solutions are generated with 
load factor N\ (see Appendix B) which leave the lower energy branch of the climb path at different 
points. The solution that just achieves v = v when E = E is chosen and then the load factor is set to 
N 2 for the transition from E to the higher energy branch. 

It is also possible to obtain an integrated solution if the small 7 assumption is not made (this was 
not possible in reference 12 because of Coriolis and Earth curvature terms). This may be of importance 
because 7 may become large in some altitude transitions. Now divide the last two of equations (17) 
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and integrate with F = 0: 


4f_ _ v(fh v -fv9) sin 7 
d'y g(N — cos 7) 


r dv _ r sn 


sm7 


cos 7 


-cfry + const. 


- In v = In (N - cos 7) + const. 
v{N — cos 7) = const. 


Applying this to both branches of the transition 

Vb(N 1 — COS 7) as V\(N\ — COS71) 
Vb(N 2 - COS 7 ) = V2(N2 - COS 72) 


or, solving for v and 7, 


V2(N 2 - cos 72) - vi(N\ - COS71) 
N 2 -N 1 


(59) 


7 B = COS 


v 2 N\ (N 2 — cos 72) — v\N 2 (Ni — cos 71) 
v 2 (N 2 - cos 72) - vi(Ni - COS71) 


Comparison of equations (58) and (59) shows that the non-small 7 solution, up, is just as easy to 
implement as the small solution, v .4. Also note that, as a check. 


lim v A = lim v B = ^1 

.Vi— +oc iVi— *00 


lim 

.Y2— *00 




lim = U2 

A^2“ + 00 


22 


NUMERICAL EXAMPLE 


The guidance algorithm developed in the previous section has been implemented in the ACSYNT 
Computer Code and used to compute near-optimal trajectories for a supersonic transport design. The 
main characteristics of the design are listed in table 1 . 

Figure 4 shows maximum thrust of the aircraft as a function of Mach number for various energy 
levels (recall that a linear throttle is not assumed), and figure 5 shows the total drag for N = 1 as a 
function of Mach for various energy levels, in the region of Mach 1. The transonic drag rise is clearly 
shown in figure 5, and this raises the possibility that there may be an instantaneous altitude transition 
in the energy climb path near Mach 1 . 

The first step of the algorithm is to find the optimum cruise point, as given by equation (42). Figure 6 
shows the optimal cruise point at each energy level throughout the flight envelope for minimum fuel 
(K\ = 0 and K 2 = 1). The optimal cruise point is interior to the flight envelope except from about 
Mach 1 .25 to Mach 1 .75, for which it is on the loft ceiling bound. 

Figure 7 shows the data of figure 6 plotted in a different way, as A x vs. Mach (see eq. (43)). This 
curve has three local minimums, each a locally optimal cruise point. One of these is a subsonic condition 
at Mach 0.95. The globally optimum point is at Mach 2.4, the highest Mach allowed. From figure 6, 
this Mach 2.4 cruise point is at an altitude of about 52,500 ft. The Mach 2.4 cruise condition has about 
a 15% higher cruise efficiency than the Mach 0.95 condition, as measured by A x ; the Mach 0.95 cruise 
point would be used for over-land flight. 

The next step in the algorithm is determining the climb path. This involves maximizing P (see 
eqs. (49) and (50)) with respect to h at energy levels from take-off to cruise. Figure 8 plots P as a 
function of Mach for various energy levels for maximum thrust and again minimum fuel. The maximum 
dynamic pressure constraint is not applied for this calculation. The value of A x used in equation (50) is 
given by equation (43) for the Mach 2.4 optimal cruise condition. It is seen that for many energy levels 
P has two or more local maxima in the vicinity of Mach 1; it is the jumping of the global maxima 
between these local maxima that causes the transonic altitude transition. 

The resulting flight path in the Mach-altitude plane is shown in figure 9. The path starts along a 
terrain limit and then climbs at almost a constant high subsonic Mach. At about 32,500 ft, it instan- 
taneously transitions to about 22,500 ft at Mach 1.25. It then continues up to the cruise point, with 
a jump to higher altitudes between Mach 1.6 and 1.8. Also shown in figure 9 is the path with the 
dynamic pressure constraint imposed. It is seen that the unconstrained path violates the constraint by 
only a small amount between Mach 1.3 and 1.7. 

Figure 10 compares the minimum fuel flight path with A x included in P in equation (50) with the 
path with the A x term omitted. The latter case corresponds to minimum fuel to climb without regard to 
a range constraint. The paths are similar except at high speed where the path with A x omitted has much 
higher dynamic pressure (there is no dynamic pressure constraint imposed). A computation was made 
to verify that including the A x term gives better performance. Referring to figure 11, the path with the 
A x term included ended with an airplane weight of 657,310 lb and a range of 836 nm. The path without 
A x ended at 690,683 lb and 349 nm. By the Brequet formula, the range covered in a cruise condition 
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is (see eq. (45)): 


■^cruise — 


v(L/D ) mp 
C mf 


At the cruise condition, v = 2323 ft/sec, ( L/D ) = 9.0, and C = 1.315 lbf ue i per hour per Thus 

for the same fuel consumed along the path with A x , the path without A x , has a range of 


349 + 


(9.0) (2323) 
(1.315) 


3600 \ 690,683 

6076/ n 657,310 


= 815 nm 


Thus the case with A x gives 21 nm. more range for the same fuel than the case without A x . 


Minimum fuel, minimum time {K\ = 0, K 2 = 1), and “minimum direct operating cost” climb 
trajectories are compared in figure 12. For the minimum cost trajectory, K\ = $500/hr and K 2 = 
$0.0626/lb; these are the values used in reference 15 for short range subsonic transports, and would 
likely need to be adjusted for supersonic long range transports. The minimum fuel and minimum time 
trajectories are quite different. The latter has no transonic altitude transition, whereas the former has 
a large one. Also, the minimum time path is much lower in altitude in the high supersonic range (the 
dynamic pressure constraint was relaxed for this calculation and would be violated by the minimum 
time path). As expected, the minimum cost path is intermediate between the other two, being more like 
the minimum time path. 

One of the principal goals of this research has been to develop an algorithm for computing the 
trajectory segments connecting the branches of the energy climb path in the transonic region, that is, 
the altitude transitions. Specifically, equation (50) was used to determine v, the value of v which is to 
be obtained when E = E (see fig. 3). An iteration is then made to determine where on the subsonic 
branch of the ECP the departure should be made to achieve this condition. The constant load factor 
integration as described in Appendix B is used to generate the flight paths. 

Figure 13 shows the transition for N\ = 0.97 and N 2 = 1.05 for the minimum fuel case in the 
altitude-Mach plane, and figure 14 shows the same path in the transonic region. The integration is 
terminated when the flight path angle is equal to the flight path angle on the supersonic branch of the 
ECP as given by equation (53). The dynamic pressure limit was ignored for this calculation. The figures 
show that there is a very close match between the altitude transition and the ECP at the termination of 
the former, and that even mild maneuvers (N\ and N 2 close to 1) give adequate transition trajectories. 

The transition trajectories for the same conditions, but using the linear estimate of v as given by 
equation (58), are shown in figure 15. Comparing figures 14 and 15 shows that the nonlinear solution 
gives a better match with the supersonic branch of the ECP than does the linear. 

The transition trajectory for a more severe load factor maneuver, JVi = 0.5 and N 2 = 1.5, is shown 
in figure 16 (these load factors would not be acceptable for a commercial transport). As compared with 
a more benign maneuver, as shown in figure 14, the transition through E occurs at a much higher 
altitude and the trajectory is much closer to E, as expected. 

Figure 17 shows the variation of energy rate, vF, as a function of Mach in the transonic region for 
the mild transition (JV X = 0.97, N 2 = 1.05). As expected, the energy rate drops when the load factor is 
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switched from 0.97 to 1.05, but never gets near zero. Also as expected, the flight path angle, 7 , at first 
decreases, and then increases when the load factor is switched as shown in figure 18; the magnitude of 
7 stays below 6 deg, making the small 7 approximation extremely good. 

The same plots are made for the more severe maneuver ( N\ = 0.5, N 2 = 1.5) in figures 19 and 
20. In this case, the energy rate becomes negative after the load factor switch and the magnitude of 
7 reaches about 22 deg, meaning that equations (59) and not equations (58) should be used for the 
calculation of the transition point. 
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CONCLUDING REMARKS 


An algorithm for optimizing supersonic transport trajectories suitable for use in an aircraft synthesis 
computer code has been developed. The algorithm has been implemented in the ACSYNT computer 
program and illustrated using a typical supersonic transport design. 

The algorithm is based on singular perturbation theory and complete time-scale decoupling of the 
energy-state version of the equations of motion (except for the fast dynamics). This results in replacing 
the functional optimization problem by a series of function optimization problems. 

The first problem is determining the optimal cruise condition. This involves a weighted sum of 
the importance of time and fuel consumption. The second problem is determining the energy-climb 
path (ECP) to the cruise condition, which involves a weighted sum of the importance of time, fuel 
consumption, and cruise efficiency. 

For the fast dynamics, a variable is introduced such that the ECP gives a consistent value of the 
flight path angle. This variable is left on the same time scale as the flight path angle and a nonoptimal 
solution of the fast dynamics using constant load factor segments is obtained. 

Numerical results for a nominal supersonic transport showed the following: (1) The optimal cruise 
point was at the highest and fastest point in the flight envelope, although there are local optimal cruise 
points at high subsonic and low supersonic speeds. (2) The ECP for the minimum fuel case had a large 
transonic altitude transition, the minimum time case had no transition, and the minimum direct operating 
cost case had a mild transition. (3) The altitude transition solutions gave good matches between the 
subsonic and supersonic branches of the ECP with operationally acceptable load factors. 

There are two obvious shortcomings of the present state of the analysis. First, the weight is 
held constant during the search for the optimal cruise point. This weight is the gross take-off weight 
according to the time-scale assumptions, but in practice could be some empirical estimate of the weight 
at the start of cruise. Second, the range during climb and descent is ignored when optimizing the cruise 
point. This is obviously a more serious problem at short ranges than for the long ranges of a supersonic 
transport. 

It is expected that both of these shortcomings could be eliminated by solving the time-scaled 
equations of motion to first order, that is, by expanding all the state variables to first order terms. This 
is the next obvious step in this research. In addition to solving these problems, the first order solutions 
will give better overall accuracy to the algorithm. 
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APPENDIX A - NOMENCLATURE 


D = drag 

E = mechanical energy 
F = normalized tangential force 
9 = gravity 
h = altitude 
H = Hamiltonian 
J = cost functional 
L = lift 
m = mass 
M = Mach number 
q = dynamic pressure 
T = thrust 
v = velocity 
x = range 

a = angle of attack 
f3 = fuel flow rate 
e = small parameter 
7 = flight path angle 
A = adjoint variable 
4> = cost function 
p = air density 



APPENDIX B - NUMERICAL INTEGRATION OF STATE EQUATIONS 


In this appendix we give the algorithms by which the state equations are integrated within ACSYNT. 
Two cases are of interest. First, the integration of the trajectory when pairs of altitude and energy 
(equivalently altitude and speed or Mach number) are given; ( Eo,ho ), (E\,hi), • • •, (Ef, hf). This is 
sometimes called path following. Second, the integration of the trajectory when the normal load factor 
N is held constant. 


Path Following 


It is assumed that all variables are known at step n — 1. These values are sought at step n, knowing 
only E n and h n . Of particular interest are the values of t n , m n , and x n . We start with equations (7), 
with 7 = 0 , written in finite difference form from step n — 1 to step n: 


Am _ 
~St ~ 


-0 


1ST = vcos7 

A E _ 

At gm 


Ah 

At 


— wsin7 



where F' = T cos a - D and, if Q is any variable. 


(Bl) 


A Q — Qn Qn — 1 



Qn ~b Qn — 1 
2 


From the first of equations (Bl), 
so that 


m n = m n _i - (3 At 
m = m n _ i - ^At 


Substituting this into the third of equations (Bl) and solving for At, 


(B2) 


At = 


Note that v, A E, h, and (3 are all known (the latter if throttle is fixed). The only quantity not known 
in equation (B3) is F' , which depends on <y. If At were known, the fourth of equations (Bl) gives 7 - 


m n -i 


vF 1 

WE 


+ 


(B3) 


7 = sin 



(B4) 
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The algorithm may now be stated as follows: 

1 . Guess a n . 

2. Compute At from equation (B3). 

3. Compute 7 from equation (B4). 

4 . Check to see if the fifth of equations (Bl) is satisfied to a suitable degree of accuracy. If not, 
select a new a n by a suitable one-dimensional search procedure and return to step ( 2 ). If satisfied, 
continue. 

5. Compute m n from equation (B2) and t n and x n from: 

t n = t n _ 1 + At 

x n = x n-l +vAtcos*f 


Phillips 35 has proposed an alternative integration scheme as follows. The third of equations (Bl) 
is now averaged directly 


' _ 2 _ f v nF n _|_ ‘^n—l^'n—1 

At 2 g y m n m n - \ 

This is then combined with equation (B2) to give 


— — 1 ) At 2 - (F n v n + F n -iv n -i + 20AE) At + 2gm n -iAE = 0 

gm n - 1 J v 7 

This is a quadratic equation to be solved for At, and replaces equation (B3) in the numerical procedure. 
As the integration step size tends to zero, these two integration schemes become equivalent. 

Constant Normal Load Factor ( N ) Paths 

In this case, AE is not a suitable integration variable because it may happen that A E < 0, which 
causes serious numerical problems. Alternative choices are At and A 7 . Because the choice Af results 
in an algorithm with three nested iterations, we follow Phillips 35 and choose A 7 . For this integration 
we do not neglect 7 . 

Because AT, A 7 , and 7 n are now known, 

~K = g(N - cos 7 ) 

is a known constant. Thus the finite difference form of the last of equations ( 8 ) is 

A7 K 
At v 
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(B5) 


Use this and equation (B2) to eliminate m n 

h n = h n —\ + 


and A t from the rest of equations (Bl). 
^ A7 sin 7 ^ + Vn-i ^ 2 


The result is 


A E = 



0 = N 


(gm n - 1 



— T sin a — X 


(B6) 

(B7) 


The algorithm is as follows: 

(1) Guess v n . 

(2) Solve for h n from equation (B5). 

(3) Compute E n = h n + and A E = E n — E n -\. 

*9 

(4) Guess a n . 

(5) Check to see if equation (B7) is satisfied. If not, select a new a n and repeat this step. If 
satisfied, continue. 

(6) Check to see if equation (B6) is satisfied. If not, select a new v n and return to step (2). If 
satisfied, continue. 

(7) Compute all other quantities of interest. 


Thus this algorithm requires a nested two parameter search, whereasjhe path following routine required 
a one parameter search. From equations (B5)-(B7) it is seen that K = 0 (N = cos 7) is not allowed. 
Should this happen, one solution is a At integration but, as mentioned earlier, this involves a three 
parameter search. 

Phillips 35 has proposed an alternative method of constant load factor integration with A7 as 
integration variable. This approach holds all variables constant at the previous step but does a second 
order integration of the altitude state equation. The increments Af and At; may be now directly 
computed from the third and fifth of equations (3): 

A( = tAt 

g(N - C0S7 n _i) 


Av = g{F n - 1 - sin7 n _i)A£ 
Differentiating the fourth of equations (3): 

h = V sin 7 + v'y cos 7 
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Using the third and fifth of equations (3) this becomes 


h = g(F n -i sin 7 „_i + iVcos 7 n _i - 1) = C 


Integrating twice: 


h = ^Cf 2 + cit 

Z 

The constant of integration ci is determined from h n _i = v n -l sin 7 n _i = c\ so that 

h n = u n _iAtsin 7 n _i + ^Ai 2 +-/Vcos 7 „_i - 1) 

with At determined as above. Reference 35 shows that this gives good numerical results. 



APPENDIX C - NECESSARY CONDITIONS FOR FAST DYNAMICS 


The state equations of the fast dynamics are the last two of equations (54): 

/ = fh v i + fv g(F - 7 ) 

7 = —(N — 1 ) 

V 

with m, x, and E (the slower states) all known constants; the control variable is a. 
From equation (23) the Hamiltonian is 


(Cl) 


H = -K x - K 2 0 + Ax v + X E vF + \ f [f h v'y + f v g(F - 7 )] 

5 

+A 7 -(iV - 1) + u i s i 

v t=l 


(C2) 


where the constraints equations ( 1 1 ) are now state constraints and must be adjoined to H with multipliers 
i/j ; the Si are assumed to be written as functions of / and E, the latter a known constant. The adjoints 
A x and Xe are known constants from the slower dynamics solutions, equations (43) and (51). From 
equation (24) the adjoint equations are 

A f = K20f - X x Vf - X E VfF - Xe vFf - Xj | (f h )f v'y + fh vp + f v gFf 
+ ifv)f g(F - 7 ) 


+ X 1^2 V A N ~ !) - A 7 i N f - E 


(C3) 


A 7 = -A f fhV + X f f v g 

where the notations of equations (5) and (40) have been used. In these equations, if Q(h,v ) is any 
function then 

Qf — Qh.hf + QvVf (C4) 


Assuming an unbounded optimal control, conditions (a) and (b) of the maximum principle give 


A £ vF a + A f f v gF a + A 7 -JV Q — 0 

-K\ -K 2 (3 + X x v + X e vF + A f [f hV1 + f v g{F - 7 )] + A 7 |(JV - 1) = 0 

From equations (2) explicit forms for F a and N a are 

F a = -E (T a cos a — T sin q — D a ) 
mg 

N a = -E (T a sin a + T cos a + L a ) 
mg 


(C5) 


(C 6 ) 
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Equations (Cl), (C3), and (C5) are used to model transitions from an initial condition to the energy 
dynamics solution (energy climb path, or ECP), from the ECP to a terminal condition, and between 
different branches of the ECP if the ECP is discontinuous. In what follows, the first case, transition 
from an initial condition to the ECP will be considered for the purpose of illustration. 

For this case, the boundary conditions on equations (Cl) and (C3) are 


m = /o 

7(0) = 70 


A/(0) = 


K\ + K2P0 - A x uo - - A 7o ^(iV 0 - 1) 


A 7 (0) = A 


70 


fho v OlQ + fv 0 9 (Fq - 70 ) 
selected to match with ECP 


(C7) 


where the second of equations (5) was used and where all quantities are known except A 70 . In summary, 
equations (Cl) and (C3) are to be integrated with control given by the first of equations (C5) subject to 
initial conditions equations (C7). 

The fast dynamics equations depend on the nature of the ECP solution because this solution 
determines the choice of variable /. If the ECP solution is an unbounded optimum, singular perturbation 
theory states that the ECP solution will be an equilibrium point of the fast dynamics, 7 ’® and the goal 
is to find a solution of the fast dynamics such that the solution approaches the ECP as t —■ ► 00 . If, on 
the other hand, the ECP is on a constraint, then the fast dynamics solution may reach the ECP in finite 
time . 36 


Some examples will now be given. If the ECP is on a terrain limit. 


/ = f(h, v) — h h m 

In this case the transformation ( h.v ) — > ( E,f ) and its inverse are given by 

E = h + -U 2 h = f 

2 9 

f = h-hm V= y/2g(E - /) 

so that 

A = 1 > fv = 0 , h f = 1 , v f = -^ 


and from equation (C4) 


Qf = Qh~ ~Qv 
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Putting these results into equations (Cl) and (C3) 

f = v 7 




~ K 2 (Ph ~ yP V ) + + X Ey F ~ ^ (7 ~ + "V 7 

2 5 

(JV -!)-*,§ (A* - f Hv) ~ E ■* (*» - %*)' 


A 7 = — XfV 

and into equations (C5) 


\gvF a + A y—N Q = 0 


v 


-Ki - K20 + X x v + A evF + \fv~f + A 7 ^(N - 1) = 0 
The initial conditions for the integration of equations (C7) are as follows: 

/( 0) = /o 

7(0) = 70 

K\ + -K 2 A) - Ax^o ~ AgupFp - A 70 g(iVo - 1) 


MO) 


u 0 70 


A 7 (0) = A 7o , selected to match with ECP 


If the ECP is on a dynamic pressure limit, the transformation ( h , v) — > (E, f) is 

E = h + -U 2 

2g 

, _ 1 2 
/ — 2 ^ — 9m 


with p = p(h ) so that 


1 2 r 

/ft = 2 PhV > /v = 


The inverse transformation is implicit. Taking differentials and using the fact that E = 


dE = dh + -dv = 0 

9 


const.: 
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Combining these equations gives 


1 o 

df = ~PhV dh + pv dv 
z 



h f= 1 


1 


§ PhV* - P9 


Then from equation (C4) 



This gives, for example. 


_ _ Ph 
Pf ~ l o 

hPh vZ ~ P9 


/. \ Phh 

W / = T-J- 


$p h v* - pg 


{f x p h v 

h f \ph,v 2 - P9 v(p- 

(f\ PH i P 

v f bhv 2 -pg v(p- 

p. | Ek 

f \ph,v 2 - pg v (p - eppj 

a _ Ph , Pv 

f \phv 2 -pg v(p- eppj 
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Equations (Cl) and (C3) become 


/ = \ph vZ 7 + P°9{P ~ 7) 


7 = ^-1) 


A f = K 2 0f - X x vj - A eV/F- X e vFj - A j [(fh)firr + -^Ph v2v fl + P v 9 F f 

5 

+ (fv) f9(F ~ 7 ) + A Avf(N - 1) - A 7 *JV, - £ Wfi 

J v v i = 1 


and equations (C5) become 


A 7 = -Xf^PhV 3 + X f pvg 


A E y Fa + A v pvgF a + A 7 -iV a = 0 


-Ph v3 l + pvg(F - 7 ) + x i~( N - 1) = 0 


— K\ — K 2 3 + X x v + A £ vF + Xj 
Initial conditions equations (C7) become 

/( 0) = fo 

7(0) = 70 

K\ ~r K 2 (3 q — X x vq — XevqFq — X 10 ^(No — 1 ) 

\ph 0 v o70 + POVog{Fo - 70 ) 

A 7 (0) = A 7o selected to match with ECP 

In all of these equations, quantities such as vj, Fj, and (fh)f are to be determined from the equations 
derived above. 

Finally, consider the case for which the ECP is an unbounded local optimum. From equation (52), 
/ in this case is 

f = Pv- -Ph 

9 
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where P is given by equation (50). Thus 


fh — Pyh g Pfih 
fv — Pvv ~ — 


Because E = const.. 


dE = dh + -dv = 0 

9 


df - ffrdh + f v dv 


so that 



1 

h-ty 


V, 


1 


fv — 



9 


Let 0 = 0(/i, v ) be defined as 


4 >- fv 


vfh 

9 


Then 


1 2“U V* 

0 = Pvv ~ “-P/i ~Phv "b ~g2^ > hh 


and equation (C4) becomes 


_ Qv- IQh 
Q f = 4, 

This explains how to compute quantities such as (3f, Ff, and (//Jy in equations (C3). Equations (Cl), 
(C3), (C5), and (C7) will not be written out explicitly for this case. 
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Table 1. Characteristics of supersonic transport 


Gross take-off weight 

753,500 lb 

Wing planform area 

5500 ft 2 

Wing span 

137.35 ft 

Leading edge sweep 

48 deg 

Aspect ratio 

3.43 

Body length 

314 ft 

Payload 

first class passengers 

30 

coach class passengers 

274 

flight crew 

2 

flight attendants 

9 

Maximum Mach number 

2.4 

Maximum dynamic pressure 

1000 psf 



Figure 1. Sketch of constraints defining flight envelope. 





Figure 2. Load factor history during altitude transition for a high performance aircraft. 



Figure 3. Sketch of an altitude transition. 
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Figure 4. Thrust vs. Mach number in the transonic region. 
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Figure 5. Drag for N 


1 vs. Mach number in the transonic region. 










MACH No. 

Figure 8. Energy rate vs. Mach number in the transonic region. 



MACH No. 

Figure 9. Energy climb path with and without dynamic pressure constraint. 
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MACH No. 

Figure 10. Energy climb path with and without A x . 



Figure 11. Sketch of trajectories with and without \ x - 
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Figure 12. Energy climb path for minimum fuel, minimum time, and minimum cost. 
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Figure 13. Transonic altitude transition for N\ = 0.97, N 2 = 1-05, nonlinear determination of v and 7. 
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MACH No. 

Figure 14. Altitude transition in the transonic region. 
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